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Abstract 



Context. The location of coronal heating in magnetic loops has been the subject of a long-lasting controversy: does it occur mostly at 
the loop footpoints, at the top, is it random, or is the average profile uniform? 

Aims. We try to address this question in model loops with MHD turbulence and a profile of density and/or magnetic field along the 
loop. 

Methods. We use the Shell-Atm MHD turbulent heating model described in Buchlin & Velli (2006), with a static mass density 
stratification obtained by the HydRad model (Bradshaw & Mason 2003). This assumes the absence of any flow or heat conduction 
subsequent to the dynamic heating. 

Results. The average profile of heating is quasi-uniform, unless there is an expansion of the flux tube (non-uniform axial magnetic 
field) or the variation of the kinetic and magnetic diffusion coefficients with temperature is taken into account: in the first case the 
heating is enhanced at footpoints, whereas in the second case it is enhanced where the dominant diffusion coefficient is enhanced. 
Conclusions. These simulations shed light on the consequences on heating profiles of the complex interactions between physical 
effects involved in a non-uniform turbulent coronal loop. 
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1. Introduction 

It is now widely accepted that magnetic loops are the basic build- 
ing block of the closed solar corona in either active regions or the 
quiet Sun, and the mechanism responsible for their heating con- 
stitutes a major unsolved problem. The energy requirements in 
the corona are quite well known, with approximately 300 WrrT 2 
needed in th e quiet Sun and 10 4 Wm ~ 2 in the active regions 
respectively (Withbroe & Noveslll977l) . It is also believed that 
photospheric motions can provide this power, with this injected 
energy being transported into the corona by waves or by the 
slower movements of the magnetic field lines. The difficulty lies 
in identifying a credible dissipation mechanism, and it is clear 
that small scales (< 1 km), such as are likely to be produced 
by turbulence, are needed. Following the work of Heyvaerts & 
Priest d 19921) and lGomez & Ferro Fontanl (1 19921) who showed 
that a self-consistent model of turbulence could account for coro- 
nal heating, many numerical simulations have been performed in 
order to study this s cenario, e.g. by 2D direct numerical simula- 
tions of MHD (e.g. iDmitruk et alJll998t[Georgoulis & Vlahosl 
119981). by cellular automata (e.g.lLu & Hamittonlll991l: Vlahos 
et al .1l995HBuchhn et alj|2003h or by shell-models (e.g. Nigro 
et al. 12004 iBuchlin & Vellill2007h . 

An important aspect of such studies which is the subject of 
this paper concerns the spatial distribution of the heating as a 
function of distance along the loop. A determination of average 
energy dissipation profile as a function of the position along a 
coronal loop would put a strong constraint on the physical mech- 
anisms involved in coronal heating. Observations hav e suggested 
that it may be uniform (as concluded for example by Pries t et ail 
1998), occur predominantly at the loop top (e.g. lRealel f2002; 



Mart ens et "ail 12002 



the footpoints (e.g. 



Schmelz & Martens! |2006l) or perhaps at 



1999; Aschwanden et al.1 
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_ lAnt iochos e t al... 
120011 : iPatsourakos et all 12004 iGudiksen & Nordlundl 12005 
However, if the heating i s due to nanoflares, simulations by 
Patsourak os & Klimchukl ((2005) show that present day obser- 
vations (including the above papers) using spectroscopic lines 
emitted around 1 MK or below can give little insight on the lo- 
calization of heating along the loop. This is because the diagnos- 
tics are measured during the cooling of the loop, long after the 
actual heating has occurred, and so depend only weakly on the 
location of the heating. The original spatial distribution of the 
heating has been smoothed out by efficient thermal conduction 
along the magnetic field. 

Furthermore, the mul ti-thermal nature of observed loops 
dSchmelz & Martensll2006l) indicates that they are composed of 
multiple unresolved sub-loops (strands), although this is contro- 
versial. In this scenario, each strand is heated intermittently and 
reaches a high temperatures (perhaps in excess of 10 7 K in an ac- 
tive regio n) before being cooled first by conduction and then by 
radiation ( jCargilll li994). As heat transport is inefficient across 
magnetic field lines (i.e. across strands), the whole process in 
each strand can usually be considered to be independent from 
the other strands. 

This paper sets out to determine the spatial location of heat- 
ing in magnetically closed structures. We simulate a thin coronal 
loop (i.e. a strand of an observable coronal loop) using a model 
of reduced-MHD (RMHD) turbulence and Alfven waves prop- 
agation. Importantly, we take into account for the first time the 
effect of a density stratification based on a one dimensional hy- 
drodynamic loop model. This represents an important first step 
in coupling energy transport along a magnetic field with MHD 
processes. The role of a non-uniform magnetic field, as well as 
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the effect of the temperature-dependence of the diffusivity coef- 
ficients are discussed. 

Our paper is organized as follows: in Sec. [2] we describe the 
models we use, in Sec. [3] we present the numerical simulations 
and their results in different cases, and in Sec.|4]we discuss the 
results further and we give some conclusions. 



2. Description of the model 

The main model used in thi s paper is the Shell-Atm model in- 
troduced in iBuchlin & VeiLI (120071) . which models MHD turbu- 
lence in a coronal loop, determining the velocity and magnetic 
field amplitudes as a function of time, position along the loop 
and wavenumber k± in the transverse direction. The shell model 
uses a mass density profile produced by an equilibrium run of 
the HydRad model which cal culates the temperature an d den- 
sity along the magnetic field (Bradshaw & Mason 2003). Both 
models are now briefly described. 



2.1. The Shell-Atm model 



The Shell-Atm codefl dBuchlin & Vellil |2007|) is a model of 
MHD turbulence in a flux tube permeated by a strong ma gnetic 
mean field B\\. It is based upon reduced-MHD (RMHD: Strauss 
Il976l) ; in this approximation of incompressible, low-p 1 , single 
fluid MHD, the magnetic field is decomposed into B — B\\ + B ± , 
where B\\ defines the z-direction (axial, or parallel direction) and 
B ± is a small perpendicular fluctuation (B x <Z By); the veloc- 
ity field is also a small perpendicular fluctuation: u — u ± with 
u± <k B\\l y/iop, where p is the mass density. 

A further simplification is made, which allows us to per- 
form simulations that are much longer and with much higher 
Reynolds numbers tha n direct numerical sim ulations of the 
RMHD equations (e.g. iDmitruk et all I2003I) : this is needed 
due to the very wide range of scales involved in coronal tur- 
bulence and because we need to perform long runs in order to 
get average profiles of the heating. The simplification consists 
of modelling the non-linear dynami cs of RMHD in planes per - 
pendicularto B\\ using a shell-model (Giuliani & Carbone 1998). 
Fourier space (which corresponds to a perpendicular section of 
the RMHD model at position z and of width t{z)) is divided into 
exponentially-spaced concentric shells k„(z) = (2n / £{z))2" , each 
containing a scalar value for the velocity and magnetic field fluc- 
tuations perpendicular to B\\ defined as u„(z) and b n (z), where 
b(z) has been normalized by ^J/j p(z). It is assumed that non- 
linear interactions are only allowed between triads (triplets) of 
neighboring shells. 

All of the shell-models distributed along B\\ are coupled 
through the Alfven waves as described by the RMHD equa- 
tions. The wave propagation takes into account the effects of a 
possible non-uniformity of the profiles of the mass density p(z) 
and of the axial magnetic fielcQ B\\(z) or of the Alfven speed 
Z?n(z) = B\\{z)l Vw(z), according to IVellil dl993l) . The model 
equations are: 

(3,±/7||5 z )Z„ ± ±iz,f5 : (lnp)±izja/7|| = -^(v t Z„ ± +v-Zj)+ 1 lJf 

(1) 



This code is publicly available from 



http : //www . arcetri . astro . it/~eric/ shell- at m/ codedoc/ 

2 Please note that giving a profile of B|[(z) is equivalent to giving a 
profile of the flux tube width £(z), due to the conservation of the mag- 
netic flux in the flux tube. 




Figure 1. Layout of the Shell-Atm model in the case of a coro- 
nal loop: shell-models are piled up along B\\ (top) and represent 
a loop, whose footpoints are anchored in the p hotosphere (bot- 



where Z^(z) = u„(z)+b„(z) are Elsasser variables, v ± = v+77, vis 
the kinematic viscosity and 77 is the magneti c diffusivity, and the 
non-li near terms are given by Eq. (4) of Giul iani & Carbona 
(119981) (with A = 2, a = 2, 6 = 5/4 and 6 m = -1/3). 

We apply this model to a coronal loop. As shown in Fig.Q] 
the perpendicular planes (each containing a shell-model) repre- 
sent cross-sections of the loop; the two end planes represent the 
photospheric footpoints of the loop; the large-scale velocity field 
Uf at the footpoints is chosen to mimic the photospheric motions: 



2mA- 



sin 2 (7rf/f) + e 2i7rB =" sin 2 (7rf/f* 



n/2)) (2) 



where z is or L (both footpoints), n corresponds the the scales 
2n/k„ of the super granulation, A an d Bare random numbers cho- 
sen as detailed in Buchlin & Vellil d2007l) . This forcing is what 
injects the energy into the model, and its amplitude is of the or- 
der of 2 km/ s. 

The model gives the energy per unit volume as Ey - 
\p 2«(I M «I 2 + I^hI 2 ) an d the heating (power of energy dissipation) 
per unit volume is calculated by computing the time derivative of 
Ey. This gives Hy = p kl(v\u„\ 2 + T]\b n \ 2 ) (where all quantities 
may depend on z): the heating per unit length is H L = H v ni 2 /4 
and the heating per unit mass is Hm = Hy jp (again, all quantities 
may depend on z). 



2.2. The HydRad model 



The HydRad code dBradshaw & Masonl2003l) is a hydrodynamic 
model of a coronal loop with a self-consistent treatment of radi- 
ation (including non-equilibrium ion populations). In this paper 
it is simply used in an equilibrium run, so as to get a hydrostatic 
profile PhydRad(z) of the mass density along the loop. We used an 
isothermal hydrostatic solution so that the density along the loop 
is stratified due to gravity alone. Since there is no thermal flux 
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(isothermal temperature profile) the background heating and ra- 
diation are in balance, and are thus identically zero. Therefore, 
we make no assumptions about the nature of the background 
heating; once coupled to the Shell-Atm model, all of the energy 
released into the loop will be from the shell model only. The loop 
length used in this run was 10 Mm. 

2.3. A simple phenomenology for the energy and heating 
profiles 

In this section we discuss what profiles of the fields, energy and 
dissipation power along the loop would be expected from a sim- 
ple phenomenology. 

A mass density profile p(z), a magnetic field profile B\\(z) and 
a profile v(z) of the dissipation coefficient (kinetic or magnetic) 
are given (in model units). The Elsasser fields are assumed to 
be both of the same order of magnitude Z (which is reasonable 
in a loop). The effects of energy dissipation and of non-linear 
interactions on the magnitudes of the fields are assumed to be 
negligible. 

The conservation of energy fluxes </>£ = +b\\E v Ji( 2 /A to- 
gether with the conservation of the magnetic flux <pg = B\\nt 2 jA 
gives then a profile Ey °c p 1 ^ 2 of the energy per unit volume. The 
dissipation power per unit volume is Hy = pvk^Z 2 oc vk\p x ^ 2 for 
a fixed k ± . However, one must bear in mind that when a wave 
travels along the loop (with non-linear interactions neglected), 
the associated k ± is affected by the variation of the width of the 
loop (in the model, the wave stays in the same shell n of the 

1 /2 

model): k ± = k n (z) = ko(z)A. n , with ko oc BJ because of the 
conservation of magnetic flux. For this reason, for a linear wave 
in the shell n, oc vB ]]P l/2 . 



3. Numerical simulations 

For all simulations, we choose model units so that the length 
of the loop is L = 1, and so as to have the below-mentioned 
values for the magnetic field Bjj and for the mass density p/ at 
the footpoints. Conversion factors to physical units for a sample 
loop are given in Table [1] but to use the same simulations for 
other loops with different physical parameters it is sufficient to 
take different values for the basic conversion factors M, L and 
T for mass, distance and time respectively. We stress that these 
conversion factors do not necessarily represent typical values of 
the actual physical quantities in the loop (the model variables 
are not necessarily of the order of unity in model units; this will 
especially be true for the diffusion coefficients). 

If PHydRad(z) is the density profile from the equilibrium 
run of the HydRad code (see Sec. 12.2b . we define R(z) = 
PhydRad(z)/PhydRad(0), with being the position of one footpoint. 
R varies between R t = 1/30 at the loop top and Rf = 1 at the 
footpoints. The profiles of mass density and Alfven speed are 
then chosen as a function of R as specified in Table [2] for each 
of the runs. The physical conditions are chosen to be the same 
at the footpoints in each run (pj = l/R t and B^ = 1), so that 
the boundary conditions of the model remain the same in the 
different runs. 

A typical value for the diffusivity coefficients v and r\ is 
5 10 -10 in model units (see Table [2] for the exact values used 
in each run), i.e. of the order of 5 10 3 m 2 /s; typical Reynolds 
and magnetic Reynolds numbers are 10 6 . This is lower than in 
the real corona (diffusivity coefficients are higher than in the real 



Table 1. Conversion factors between model units and S.I. phys- 
ical units: expression in the general case and value for a typical 
modelled loop. 



Conversion factor to S.I. physical units 



Quantity 


Expression 


Typical value 


Mass 


M 


10" kg 


Length, distance 


L 


10 7 m 


Time 


T 


10 s 


Velocity, Alfven speed 


LT-' 


10 b ms-' 


Magnetic field 


^Afi/2 L -i/2 r -i 


1.1 10- 3 T 


Mass density 


ML' 3 


10-' 2 kgm~ 3 


Diffusivity 


L 2 T- [ 


10 13 m 2 s-> 


Energy per unit volume 


ML- [ T- 2 


1 Jm 3 


Power per unit volume 


ML-'T- 1 


10 1 Wnr 3 



Table 2. Summary of the parameters of the modelled loop for 
different runs, in model units. 



Run 


pR, 




All 


v/lO" 10 


;y/10- 10 


(u) 


1 


l 


1 


5 


5 


(a) 


R 




1 


5 


5 


(b) 


R m 


g-W 


1 


5 


5 


(c) 


R l/4 


R -m 


1 


5 


5 


(A) 


R 


R -m 


R m 


5 


5 


(B) 


R l/2 




R m 


5 


5 


(C) 


Rl/4 


U-l/16 


Rl /U 


5 


5 


(P) 


R 


R -l/2 




0.5(R/R,r V2 


0.5(R/R,) 3/2 


(q) 


R 


R -U2 




50(R/R,r V2 


0.5(R/R,) 3/2 


(r) 


R 


R -U2 




W 2 (RIR,)- V2 


1Q- 2 (R/R,) 3 > 2 


(s) 


R 


R -U2 




W 2 (RIR,)- 112 


10-\R/R,) 312 


(t) 


R 


R -l/2 




io 2 (tf/.R,r 7/2 


10- 6 (R/R,) 3/2 



corona) but it represents a huge improvement over direct numer- 
ical simulations of MHD. 

We use the model with 1000 planes piled up along B\\, and 
each of these has 18 shells (n = 0, . . . 17). Each of the simula- 
tions is run during 500 units of time (after an initial phase of en- 
ergy growth), which corresponds to the order of 100 large-eddy 
turn-over times. Profiles of the energy and of the dissipation per 
unit volume are averaged over the duration of the runs. 

3. 1 . Test case: uniform loop 

We present run (u) as a test case, where the mass density, mag- 
netic field, Alfven speed and diffus ion coefficients are all uni- 
form (as in iBuchlin & Velhl 120071) . Figure [2] shows that the 
time-average profile of the total energy per unit volume is uni- 
form along the loop, with the magnetic energy dominating. The 
kinetic energy density is also uniform except near the foot- 
points, where it gets lower because of the boundary conditions 
(uf <k 

The energy dissipation (right panel of Fig. [2j is slightly non- 
uniform: although the phenomenology of Sec. 12.31 suggests a 
uniform profile, there is a drop near the footpoints and an en- 
hancement near the loop top; the former can be attributed mainly 
to the contribution of the profile of kinetic energy dissipation, 
and the latter to the profile of the magnetic energy dissipation. 
The ratio of total dissipation power in the central 10% of the 
loop length to the total dissipation power in both 5% of the loop 
length at the footpoints is 1.5. 

A second observation is that the ratio of magnetic to kinetic 
energy dissipation is lower than the ratio of the magnetic to ki- 
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Figure 2. Average profiles of the energy (left) and dissipation 
power (right) per unit volume, as a function of the position z 
along the loop, for run (u). The dotted lines are the kinetic en- 
ergy or dissipation power, the dashed line the magnetic energy 
or dissipation power, and the plain line is the total. 



T 




2 3 4 5 6 

log 10 kp er p 

Figure 3. Time-averaged kinetic (dots) and magnetic (dashes) 
perpendicular ID energy spectra in the loop for run (u). A line 
of slope -5/3 is shown for reference. These spectra are represen- 
tative of the spectra obtained for the other runs, although there 
can be some variations between runs. 

netic energy. This is because most magnetic energy is concen- 
trated at large scaled as can be seen on the spectra of Fig. [3] 
which do not dissipate as much as the smaller scales. Then the 
relative weight of the kinetic dissipation in the total dissipation 
is increased compared to the weight of the kinetic energy in the 
total energy; as the kinetic dissipation profile is lower at the foot- 
points (because of the boundary condition with Uf <K B\\), this 
explains the non-uniform shape of the total dissipation profile. 
This kind of effect needs to be kept in mind when analyzing the 
other runs. 

3.2. Loop with a non-uniform mass density 

We now consider a more realistic loop, with a non-uniform mass 
density determined by the HydRad model. The ability to take 
into account a non-uniform mass density (as well as a variable 
axial magnetic field: see later) is an important advan tage of the 
Shell- Atm mode l when compared to the model of iNigro et al.l 
(2004) as used bv lReale et al.l d2005h . Run (a) is performed using 
a mass density profile p(z) = R(z)/R t (Fig. In runs (b) and (c) 

3 This is a well-known effect coming f rom the resonances in the so- 
lution of the linearized RMHD equations ( Milano et al. 1993). 
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0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 

z z z 

Figure 4. Profiles of mass density p, Alfven speed b\\ and gra- 
dient of Alfven speed &b\\l dz for run (a), as a function of the 
position z along the loop. 

we use R l ^ 2 /R t and R l ^ 4 /R, for the mass density profiles, so as 
to assess the effect of smoother profiles. 

Loops seen with Yohkoh and TRACE have only a slight 
width enhancement at the top relative to the footpoints 
dKlimchuk et al.l Il992t iKlimchukl l2000l ; IWatko & Klimchukl 
2000). It is thus reasonable to assume a uniform axial mag- 
netic field along such loop. In the model we set By = b\\ yfp~, 
with b\\ = 1 (jUo = 1 in model units). The Alfven speed is 
then b\\{z) = b\'(p{z)lp,) 1 ^ 2 , a function of z, as shown in Fig. [4] 
Consequently the waves are partially reflected due to the Alfven 
speed gradient . 

The resulting time-averaged profiles of the energy per unit 
volume Ev(z) and of the heating per unit volume Hv(z) are 
shown in Fig. [5] The striking result is that in all these cases, in- 
cluding run (a) which has steep gradients of density and Alfven 
speed, the effect of the non-uniform density on the heating pro- 
files is very limited. For run (a) for example, while the phe- 
nomenology of Sec. |2.3| predicts a dependence in p 1 ^ 2 for both 
energy density and dissipation, i.e. a ratio « 5.5 between the 
footpoints and the loop top, the actual ratio is only about 1.3 
between the highest values (near the footpoints) and the lowest 
values, for both energy and and dissipation. For runs (b) and (c) 
(which have shallower density gradients), the profiles are even 
closer to being uniform. 

We performed fits of each of these profiles to a power-law 
p a of the mass density (Table [3j. The domain for the fit was 
restricted to the central 90% of the loop, to avoid effects due 
to the boundary conditions. In cases (a) and (b), the fit of the 
heating profile gives an index an close to for the power-law 
p a " (taking into account the error bars of the fit; see Table [3]for 
the detailed results of the fits). Each of the individual profiles 
cannot be considered as uniform (as would be the case if the fit 
was perfect and with an - 0), but they are still considerably 
more uniform than that expected from the steep gradient of the 
mass density p. 

The fit of the energy profiles gives an index 7.1 10 -2 < as < 
9.4 10~ 2 (in all three cases, taking into account the error bars). 
The profile of energy is thus also much more uniform than ex- 
pected from the phenomenology (a = 1 /2). This is surprising as 
it means that, when going from the footpoints to the loop top, the 
energy decreases slower than what would be expected from the 
conservation of the energy flux (under the assumptions detailed 
in Sec. I2.31 >. The discrepancies between the phenomenology and 
the results are discussed further in Sec.|4] 

Although the comparison of the heating in loops of different 
lengths is out of the scope of this paper, we have also performed 
a simulation of a 40 Mm-long loop. In this case a uniform chro- 
mospheric temperature yields unrealistically low coronal densi- 
ties, so the loop has been divided into a 20 000 K chromosphere 
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Figure 5. Average profiles of the energy (left) and dissipation 
power (right) per unit volume, as a function of the position z 
along the loop, for runs (a), (b) and (c) (from top to bottom). 
The dotted lines are the kinetic energy or dissipation power, the 
dashed line the magnetic energy or dissipation power, and the 
plain line is the total (same as in Fig. |2j». 



Figure 6. Average profiles of the energy (left) and dissipation 
power (right) per unit volume, as a function of the position z 
along the loop, for runs (A), (B) and (C) (from top to bottom). 
The dotted lines are the kinetic energy or dissipation power, the 
dashed line the magnetic energy or dissipation power, and the 
plain line is the total (same as in Fig.|2]l. 



Table 3. Summary of the indices cue and an for the power-law 
fits Ey oc p" E and Hy p a ", (1) expected from the phenomenol- 
ogy of Sec. |2.3| and (2) obtained from the numerical simulations, 
for the runs with uniform dissipation coefficients. The fits from 
the simulations are performed on time-averaged profiles and on 
the central 90% of the loop. 



Run a E (1) a E (2) a H (1) a H (2) 



(a) 


1/2 (8.50 


+ 0.05)10- 2 


1/2 


(7.39 + 


0.09) 10- 2 


(b) 


(7.14 


+ 0.02)10~ 2 




(-8.98 


t0.15)10" 2 


(c) 


(9.42 


+ 0.01)10- 2 




(-2.85 


t 0.03)10-' 


(A) 


(3.15 


+ 0.00)10-' 


3/4 


(4.59 + 


0.01)10-' 


(B) 


(3.18 


+ 0.00)10"' 




(3.81 + 


0.01)10-' 


(C) 


(3.05 


±0.00)10-' 




(1.96 + 


0.04)10-' 



and, from the points where the density has been divided by 100 
compared to the footpoints, a 1 MK corona. Although the ratio 
of footpoint to loop top densities is higher than in the other sim- 
ulations presented in this paper, the average heating profile re- 
mains quasi-uniform (an = (3.13 + 0.66)10- 2 ). Therefore, when 
the heating model Shell-Atm is considered alone in a hydro- 
static loop (i.e. neglecting the feedback effects of the cooling 
processes), longer loops do not seem to behave differently than 
shorter loops. 



3.3. Loop with non-uniform mass density and magnetic field 

In very large structures like coronal streamers, it is likely that the 
magnetic field gets weaker with altitude. We model this variation 
through a dependence on the mass density: the axial magnetic 
field of the loop is chosen to be B n (z) = b\\ yj '/j p(z)(p(z) IptT 1 14 
(the exponent 1/4 has been chosen so as to minimize both the 
gradients of magnetic field and of Alfven speed for a given den- 
sity profile). The Alfven speed is then b\\(p(z)/pt) 1 ^. 

Runs (A), (B) and (C) are performed with this magnetic field 
profile, and the same density profiles as runs (a), (b) and (c) re- 
spectively (see Tabled. For run (A) for example, the ratio of the 
magnetic field between the footpoints and the loop top is then 
2.3. 

The resulting time-averaged profiles of the energy and the 
heating per unit volume are shown in Fig. [6] In all cases (and 
especially (A), which has the steepest density and magnetic field 
gradients) the heating is enhanced at the footpoints. The results 
of the power-law fits of the heating profiles as a function of the 
mass density are also shown in Tabled we get 0.19 < &h < 0.46 
in all three cases (taking into account the error bars) whereas the 
phenomenology of Sec. |2.3| predicts cr# = 3/4. 

This enhancement of the dissipation power at the footpoints 
can be seen as a consequence of the higher energy per unit vol- 
ume there. Indeed, the power-law fits of the energy profiles give 
0.30 < cue < 0.32 in all three cases, as seen in Table [3] (taking 
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Figure 7. Profiles of the dissipation coefficients v (kinematic vis- 
cosity; dots) and 77 (magnetic diffusivity; dashs) used in runs (p), 
(q), (r), (s) and (t) (from left to right and top to bottom). 



into account the error bars), whereas the phenomenology pre- 
dicts oe = 1/2. Again, there is more energy at the loop top 
than what would be expected from the conservation of the en- 
ergy flux. 



3.4. Loop with non-uniform mass density and diffusivity 
coefficients 

Up to now the diffusion coefficients have been assumed to be 
uniform, but this is not the case if the temperature profile is not 
uniform. So as to take this effect into account, we consider a 
different situation than in the previous sections, with the same 
density profile p(z) = R(z)/R t as runs (a) and (A), but assuming 
now that the loop has a uniform pressure and that the law of per- 
fect gases holds, hence a variation of a factor Rf/Rt = 30 of the 
temperature between the footpoints and the loop top. With these 
assumptions, the variations of the diff usivity coe f ficient s in the 
magnetized plasma as obtained from Braginskij| dl965h are in 
T 1 ^ 2 for the kinematic viscosity and in r~ 3 ^ 2 for the magnetic 
diffusivity, where T is the plasma temperature (which is only 
known in our model from the density p and the assumptions we 
make in this section). We thus choose profiles v = v,R~ 112 and 
77 = rj r R 3 ^ 2 for the kinematic viscosity and for the magnetic dif- 
fusivity respectively, where v, and 77, are the values of v and 77 
that we have chosen at the loop top. 

We perform runs (p), (q), (r), (s) and (t) at different magnetic 
Prandtl numbers at loop top: Pr' M = v t jt] t = 1, 10 2 , 10 4 , 10 6 and 
10 8 respectively (see Table|2|. With the profiles we use, the mag- 
netic Prandtl number is (Rf/R t ) 5 times lower at the footpoints 

than at the loop top, i.e. 2.4 10 7 times lower: it is PrL = 4 10~ 8 , 
4 1(T 6 , 4 1(T 4 , 4 1(T 2 and 4 for runs (p), (q), (r), (s) and (t) re- 
spectively (see Fig. |7). These ranges of magnetic Prandtl num- 
bers and of diffusivity coefficients, which are impossible to ex- 
plore with direct numerical simulations, are easily managed by 
the Shell-Atm model. 

Compared to run (a), which has otherwise the same param- 
eters, the profiles of the kinetic and magnetic energies, seen in 
Fig. [8] are barely affected: such small dissipation coefficients, 
which allow for a wide turbulent inertial range to develop, have 
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Figure 8. Average profiles of the energy (left) and dissipation 
power (right) per unit volume, as a function of the position z 
along the loop, for runs (p), (q), (r), (s) and (t) (from top to 
bottom). The dotted lines are the kinetic energy or dissipation 
power, the dashed line the magnetic energy or dissipation power, 
and the plain line is the total (same as in Fig. |2j. 



almost no influence on the energy content as a function of the 
position along the loop, and the energy does not decrease signif- 
icantly with altitude as a result of the dissipation. 

However, as expected, the profiles of both the kinematic vis- 
cosity v and the magnetic diffusivity 77 have an effect on the ki- 
netic {HV) and magnetic (H v ) dissipation powers (right panels 



of Fig. |8): the profile of H' v is mainly affected by the profile of v 
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and the profile of H b y is mainly affected by the profile of r\. The 
combination of both gives a heating profile which is enhanced 
either at the footpoints (runs (p) and (q)), at the loop top (runs 
(s) and (t)) or both (run (r)). Heating enhancements at the foot- 
points are mainly due to magnetic energy dissipation, while a 
heating enhancement at the loop top is mainly due to kinetic en- 
ergy dissipation. Note also that, as the kinetic energy is smaller 
in all cases than the magnetic energy, the ratio of the kinetic over 
magnetic dissipation is lower than what would be expected from 
the mere ratio of the kinematic viscosity over the magnetic dif- 
fusivity. 

The fits of the kinetic (resp. magnetic) dissipation profiles to 
a power-law p a are in most cases broadly consistent with the ex- 
pected value a — -7/2 (resp. 3/2), although there are important 
differences between the different runs. The ratios between the ki- 
netic (resp. magnetic) dissipation power near the footpoints and 
at the loop top are also broadly consistent with the same ratio for 
v (resp. iff. Given the parameters of the model, this ratio can be- 
come very large but it is again easily managed by the Shell- Atm 
model. 

4. Discussion 

We have presented a set of simulations of a coronal magnetic 
loop containing a highly turbulent MHD plasma, with Reynolds 
numbers up to 10 6 , a variation between the loop top and the foot- 
points by a factor of up to 30 for the density and up to the order 
of 10 7 for the magnetic Prandtl number. These ranges of physical 
parameters cannot be reached by direct numerical simulations. 

The spatial distributions of the energy and of the heating as a 
function of the position along the loop have been obtained. The 
profiles of energy do not have a strong dependence on the posi- 
tion along the loop (runs (a) to (c)), except when the flux tube 
is expanding (runs (A) to (C)): in this case the energy per unit 
volume is higher near the footpoints, where the magnetic field is 
higher. The profiles of the heating (the dissipation of energy per 
unit volume) also have a low dependence on position when only 
the density varies. If the flux tube is expanding, the behavior of 
the heating profile is quite similar to that of the energy profile: 
the heating is higher near the footpoints. 

When we look at runs (a) to (C) the profiles of energy and 
heating are shallower than what is expected from the simple lin- 
ear phenomenology introduced in Sec. 12.31 We now try to ex- 
plain the origin of these discrepancies. First, one might expect 
that the dissipation occurring during the wave propagation to the 
loop top would make the profile of energy Ey deeper instead 
of shallower; however, the comparison of the energy profiles in 
runs (p) to (t) (Fig. [8] left) with run (a) shows that this effect 
remains unnoticeable. The origin of the energy profile is still un- 
clear. 

Assuming now an energy profile Ey oc p" E (with oe coming 
from the fits of the results of the different runs), the dissipation 
profile in a given shell n is Hy oc B\\p aE , i.e. an ~ 0.1 for 
runs (a) to (c) and an w 0.55 for runs (A) to (Cfl This is al- 
ready closer to the simulations results (Table [3) than the purely 
phenomenological predictions of 1/2 and 3/4 respectively. The 
remainder of the difference could be explained by the effect 
observed in run (u): the ratio of the kinetic dissipation profile 
(which is low at footpoints due to the velocity boundary con- 
dition at the photosphere) over the magnetic dissipation profile 

4 It can also be noted that these new exponents are close to the ex- 
ponents a H = a nd ag = 1/2 that would be ob tained from a heating 
proportional to Bl (Gudiksen & Nordlund 2005a). 



is increased compared to the ratio of the kinetic energy profile 
over the magnetic energy profile because of the concentration of 
magnetic energy at the larger scales in the loop. Thus the slow 
motions of the loops footpoints anchored in the photosphere and 
the build-up of magnetic energy mainly at large scales have a 
decisive influence on the determination of the profiles of heating 
in coronal loops. 

Another, more subtle, effect is related to the profiles of the 
nonlinear and dissipation timescales. For a given k ± , the non- 
linear time scale is tnl = l/kj_Z and the dissipation time scale 
is t v = \jvk\. Taking into account the energy profile and the 
variation of k ± for a given shell n (like in Sec. 12. 3K the non- 
linear time scale is r^ L oc p (1_Q,£)/2 /B^ 2 and the dissipation time 

scale is Ty oc l/vB\\ for a given shell n. This means that for 
runs (a) to (c), the dissipation time scale does not depend on po- 
sition while the nonlinear timescale (« p 045 ) is shorter at the 
loop top: the transfer of energy to the small scales is more effi- 
cient, and as a consequence, the dissipation power is higher at 
the loop top than what is expected from the phenomenology of 
Sec. 12.31 This is in accordance to what is seen in the simulation 
results (a shallower profile of dissipation). For runs (A) to (C) 
the nonlinear (oc p 023 ) and dissipation (oc p~ I/4 ) timescales vary 
in opposite directions: the nonlinear terms are more efficient at 
the loop top, but they need to bring the energy further in the spec- 
trum. However, the different ^-dependence of both competing 
timescales means that the enhanced efficiency of the nonlinear 
transfer "wins" over the displacement of the dissipation scale: 
again, as seen in the simulations, the heating profile is shallower 
than expected before. 

Runs (p) to (t) point out the additional role of the variation of 
the diffusion coefficients as a function of temperature (and thus 
as a function of position along the loop), which seems to have 
been overlooked in previous works. The profile of the kinetic 
(magnetic) energy dissipation follows approximately the pro- 
file of the kinetic (magnetic) diffusion coefficient respectively. 
The diffusion coefficients have thus a direct influence on the en- 
ergy dissipation, and as they strongly vary with temperature, this 
yields a strong variation of the dissipation power along the loop. 
Furthermore, as both diffusion coefficients vary in opposite ways 
with temperature, the heating can be enhanced at the footpoints, 
at the loop top or both, depending on the magnetic Prandtl num- 
bers in the loop. 

Although these results seem to be straight forward, this had 
never been modelled in the context of solar coronal heating, 
as the precise values of v and 77 are usually considered unim- 
portant to get averages of the heating: the main argument for 
that is that energy is dissipated at the end of the inertial range, 
at whichever scale (or wavenumber) this end is (Galsgaard & 
Nordlund 119961) . However, these values actually have an im- 
portance because (1) they discriminate between low- and high- 
Reynolds-number physics, with effects like intermittency ap- 
pearing only at very high Reynolds numbers; (2) when waves 
propagate between regions with different diffusion coefficients, 
we have shown in this paper that the profiles of these coefficients 
contribute to the determination of the profile of the heating. 

These simple simulations shed light on the consequences on 
heating profiles of the complex interactions between some of 
the physical effects that come into play in a non-uniform turbu- 
lent coronal loop: conservation of energy and magnetic fluxes, 
accumulation of large-scale magnetic field in a loop submitted 
to photospheric motions, nonlinear effects, and non-uniformity 
of the viscosity and magnetic diffusivity coefficients. However, 
many processes have not been included in these simulations, 
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such as heat transport (conduction and radiation), flows, grav- 
ity, MHD processes overlooked by the shell-models, and kinetic 
effects. Future developments may include the heat transport and 
the subsequent computation of emitted radiation, especially in 
UV. In addition to completing the picture of the nanoflare-like 
processes involved in heating and cooling of loops, this will al- 
low the comparison of such models to data from new instru- 
ments (such as Hinode/EIS) that observe the plasma at higher 
temperatures than previous instruments, at earlier stages of heat- 
ing events. 
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